#include <iostream>
#include <string>
#include <vector>
#include <cstdlib>

#include "filter.h"

#include "Image.H"
#include "IO.H"
#include "Util.H"

using namespace pip;
using namespace std;

/// MEDIAN FILTERING \\\

/**
 * Performs median filtering on an image, using mirrored borders.
 *
 * @param image Image to filter
 * @param n 
 * @param squareMask
 */
void medianFilter(Image<unsigned char> &image, const int &n, const bool &squareMask = true)
{
	//check if kernel is odd
	if (n%2 == 0){
		cerr << "given width is not odd. stopping medianFilter..." 
			 << endl;
		return;
	}
	
	int range((n-1)/2);
	Image<unsigned char>  temp  = image.bordered(range,border_mirror);
	vector<unsigned char> values;
	//vector<unsigned char>::iterator mid;
	int median = ( squareMask ? (n*n/2+1) : n); // median's place in ordered sequence

	//for all  pixel in image
	for (int x=0; x < image.dim(0); ++x){
		for (int y=0; y < image.dim(1); ++y){
			//fill values within kernel into vector
			if(squareMask){
				for(int dx = -range; dx <= range ; ++dx){
					for(int dy = -range; dy <= range ; ++dy){
						values.push_back(temp(x+dx,y+dy));
					}
				}
			}
			else{
				for(int dx = -range; dx <= range ; ++dx){
					values.push_back(temp(x+dx,y));
				}
				for(int dy = -range; dy <= range ; ++dy){
					if(dy!=0) values.push_back(temp(x,y+dy)); //middle element is already in vector
				}	
			}
			nth_element(values.begin(),values.begin()+ median,values.end());
			image(x,y) = values[median];
			values.clear();
		}
	}


}

/// SHARPENING \\\

/**
 * Initializes the Laplacian kernel
 *
 * @param kernel The Laplacian kernel
 * @param scale
 */
void initLaplacianKernel(Image<float> &kernel, const float &scale)
{
	kernel(0,0)= 1*scale; kernel(0,1)= 1*scale;  kernel(0,2)= 1*scale;
	kernel(1,0)= 1*scale; kernel(1,1)= -8*scale; kernel(1,2)= 1*scale;
	kernel(2,0)= 1*scale; kernel(2,1)= 1*scale;  kernel(2,2)= 1*scale;
}

/**
 * Sharpens the image by subtracting a Laplacian with factor C.
 * See also 5.35 in Sonka
 *
 * @param image Image to sharpen
 * @param C Factor
 */
void sharpen(Image<float> &image, const float &C)
{
	Image<float>  temp = image;
	Image<float>  temp2 = image;
	Image<float>  mask(Dimension(3,3));
	initLaplacianKernel(mask,C);

	filter(temp, temp2, mask, border_mirror);

	for(int x=0; x < image.dim(0); ++x){
		for(int y=0; y < image.dim(1); ++y){
			image(x,y) = temp(x,y) - temp2(x,y);
		}
	}
}

/// UNSHARP MASKING \\\

/**
 * Performs unsharp masking on an image.
 * 
 * @param image Image to filter
 * @param sigma
 * @param a
 */
void unsharpMasking(Image<float> &image, const float &sigma, const float &a)
{
	Image<float>  blurred  = image;
	Image<float>  highpass = image;

	gauss(image, blurred, sigma);

	for(int x=0; x < image.dim(0); ++x){
		for(int y=0; y < image.dim(1); ++y){
			highpass(x,y)=image(x,y)-blurred(x,y);
		}
	}

	for(int x=0; x < image.dim(0); ++x){
		for(int y=0; y < image.dim(1); ++y){
			image(x,y) = blurred(x,y) + highpass(x,y)*a;
		}
	}
}

/// DENOISE \\\

/**
 * Generates a random number in the range [min, max]
 *
 * @param min Lower range value
 * @param max Upper range value
 * @return Random number
 */
float randomNumber(const float &min, const float &max)
{
	if(max>=min){
	
		return (rand()/RAND_MAX * (max-min) +min); 
	}
  return 0.0f;
}

/**
 * Generate a random number from a Gaussian distribution N(mean, sigma)
 *
 * @param mean
 * @param sigma
 * @return Random number
 */
float randomGaussian(const float &mean, const float &sigma)
{
  float U1, U2, V1, V2, S;

  do {
    U1 = randomNumber(0, 1);
    U2 = randomNumber(0, 1);
    V1 = 2 * U1 - 1;
    V2 = 2 * U2 - 1;
    S = V1 * V1 + V2 * V2;
  } while (S >= 1);

  float X = sqrt(-2 * log(S) / S) * V1;

  return (mean + sigma * X);
}

/**
 * Degrade an image with additive Gaussian noise.
 *
 * @param image Image to degrade
 * @param sigma
 */
void degradeGauss(Image<unsigned char> &image, const float &sigma)
{
  for (int y = 0; y < image.dim(1); ++y) {
    for (int x = 0; x < image.dim(0); ++x) {
      float randGauss = randomGaussian(0, sigma);

      int sign = randGauss < 0 ? -1 : 1,
        greyVal = (int)floor(abs(randGauss) + 0.5),
        noisyValue = image(x, y) + sign * greyVal;

      if (noisyValue < 0) {
        image(x, y) = 0;
      } else if (noisyValue > 255) {
        image(x, y) = 255;
      } else {
        image(x, y) = (unsigned char) noisyValue;
      }
    }
  }
}

/**
 * Degrade an image with impulse noise. Fraction denotes the
 * fraction of pixels affected by corruption.
 *
 * @param image Image to degrade
 * @param fraction
 */
void degradeImpulse(Image<unsigned char> &image, const float &fraction)
{
  for (int y = 0; y < image.dim(1); ++y) {
    for (int x = 0; x < image.dim(0); ++x) {
      float randUniform = randomNumber(-1, 1);

      if (abs(randUniform) < fraction) {
        image(x,y) = randUniform < 0 ? 0 : 255;
      }
    }
  }
}

/**
 * Average squared difference ...
 *
 * @param input ...
 * @param secondInput ...
 * @return ...
 */
float averageSquaredDifference(const Image<unsigned char> &input, const Image<unsigned char> &secondInput)
{
	if(input.dim(0) == secondInput.dim(0) && input.dim(1) == secondInput.dim(1)){
		double temp;
		for (int y = 0; y < input.dim(1); ++y) {
			for (int x = 0; x < input.dim(0); ++x) {
				temp = (double)(input-secondInput)*(input-secondInput);
			}
		}
		return (temp/input.dim(0)*input.dim(1));
	}
	cerr << " Images of different size can not be compared " 
		 << endl;
	return 0.0f; //
}

/**
 * Degrades an image with noise (additive Gaussian or Impulse noise),
 * and restores it again (by Gaussian blurring or median filtering).
 *
 * @param image Input image
 * @param degradedImage Output degraded image
 * @param degradedFilteredImage The 'restored' degraded output image
 * @param s1 Parameter for degrading (sigma for Gaussian noise or fraction for impulse noise)
 * @param s2 Parameter for restoring (sigma for blurring or n for median filtering)
 * @param filterMethod True for Gaussian blurring, false for median filtering
 */
void denoise(const Image<unsigned char> &image, 
             Image<unsigned char> &degradedImage, 
             Image<unsigned char> &degradedFilteredImage, 
             const float &s1, 
             const float &s2, 
             const bool &filterMethod,
             void (*additiveNoise)(Image<unsigned char>&, const float &))
{
	degradedImage		  = image;
	additiveNoise(degradedImage,s1);
	degradedFilteredImage = degradedImage;
	if(filterMethod){
		gauss(degradedImage, degradedFilteredImage, s2);
	}else{
		medianFilter(degradedFilteredImage, s2);
	}
}

/// CANNY EDGE DETECTION \\\

/**
 * Gradient magnitude ...
 *
 * @param in Input image
 * @param dx Output of X-derivative
 * @param dy Output of Y-derivative
 * @param gm Output of Gradient Magnitude
 * @param sigma
 */
void gradientMagnitude(const Image<unsigned char> &in, Image<float> &dx, Image<float> &dy, Image<float> &gm, const float &sigma)
{
  // Add your code here!
}

/**
 * Canny edge detector ...
 *
 * @param in Input image
 * @param out Output of gradient magnitude for pixels that survive non-maximal suppression
 * @param sigma Scale for the gradient magnitude
 */
void cannyEdgeDetector(const Image<unsigned char> &in, Image<float> &out, const float &sigma)
{
  // Add your code here!
}

/**
 * Number of edge voxels ...
 *
 * @param filteredImage ...
 * @return ...
 */
int nrOfEdgeVoxels(const Image<float> &filteredImage)
{
  // Add your code here!

  return 0;
}


/// BODY OF THE PROGRAM \\\
/// You do not need to edit anything below this line!!
int applyMedian (int argc, char **argv)
{
	if (argc < 8) {
		error("", "Usage: MoreFiltering -operation Median -i <input> -o <output> "
		"-width N (and optional) -kernel S(quare) or P(lus)");
	}

	// Store commandline args in stl container
	vector<string> args(argv, argv + argc);

	// Variables to be read from the commandline
	string program, inputfile, outputfile, kernel("S");
	int width(3);

	// Read the arguments.
	getProgramName(args, program);
	getScalarArg(args, "-width", width);
	getScalarArg(args, "-i", inputfile);
	getScalarArg(args, "-o", outputfile);
	getScalarArg(args, "-kernel", kernel);

	bool square(true);
  if(kernel == "P") {
	   square = false;
  }

	cout << "Applying Median Filtering on " << inputfile << " using " << (square ? "square" : "plus-shaped") << " kernel with size " << width << "." <<  endl;
	Image<unsigned char> image;
	if(!importFile(image, inputfile)) {
		cerr << "Error loading " << inputfile << endl;
		return 1;
	}

	//getting min and max for assignment
	unsigned char min, max;
	getMinMaxImage(image, min, max);
	cout << "min of image before filtering: "<< (int) min << endl;
	cout << "max of image before filtering: "<< (int) max << endl;

	medianFilter(image, width, square);

	getMinMaxImage(image, min, max);
	cout << "min of image after filtering: "<< (int) min << endl;
	cout << "max of image after filtering: "<< (int) max << endl;

	if(!exportFile(image, outputfile)) {
		cerr << "Error writing " << outputfile << endl;
		return 1;
	}

	cout << "Successfully saved output as " << outputfile << "." << endl;

	return 0;
}

int applyUnsharp(int argc, char **argv)
{
	if (argc < 5) {
		error("", "Usage: MoreFiltering -operation UnsharpMask -i <input> -o <output> [-sigma s] [-A a]");
	}

	// Store commandline args in stl container
	vector<string> args(argv, argv + argc);

	// Variables to be read from the commandline
	string inputfile, outputfile, program;
	float sigma(2),
    a(1.5);

	// Read the arguments.
	getProgramName(args, program);
  getScalarArg(args, "-i", inputfile);
  getScalarArg(args, "-o", outputfile);
  getScalarArg(args, "-sigma", sigma);
  getScalarArg(args, "-A", a);

  cout << "Applying Unsharp Masking on " << inputfile << " using sigma: " << sigma << " and A: " << a << "." <<  endl;

  //Image<float> image;
  Image<unsigned char> imageUC;
  if (!importFile(imageUC, inputfile)) {
    cerr << "Error loading " << inputfile << "!" << endl;
    error(program, "Error reading image");
  }

  Image<float> image(imageUC.dimension());
	convertUnsignedCharToFloat(imageUC, image);

	//getting min and max for assignment
	float min, max;
	getMinMaxImage(image, min, max);
	cout << "min of image before filtering: " << min << endl;
	cout << "max of image before filtering: " << max << endl;

	unsharpMasking(image, sigma, a);

	getMinMaxImage(image, min, max);
	cout << "min of image after filtering: " << min << endl;
	cout << "max of image after filtering: " << max << endl;

  convertFloatToUnsignedChar(image, imageUC);

  if (!exportFile(imageUC, outputfile)) {
    error(program, "Error writing file");
  }

  cout << "Successfully saved output as " << outputfile << "." << endl;
  return 0;
}


int applySharpening(int argc, char **argv)
{
  if (argc < 6) {
    error("", "Usage: MoreFiltering -operation Sharpening -i <input> -o <output> [-factor F]");
  }

  // Store commandline args in stl container
  vector<string> args(argv, argv+argc);

  // Variables to be read from the commandline
  string inputfile, outputfile;
  string program;
  float C(0.7f);

	// Read the arguments.
  getProgramName(args, program);
  getScalarArg(args, "-i", inputfile);
  getScalarArg(args, "-o", outputfile);
  getScalarArg(args, "-factor", C);

  cout << "Applying Sharpening on " << inputfile << " using C = " << C << "." <<  endl;
  Image<unsigned char> imageUC;
  if (!importFile(imageUC, inputfile)) {
    cerr << "Error loading " << inputfile << "!" << endl;
    error(program, "Error reading image");
  }

	Image<float> image;
  image.resize(imageUC.dimension());

  convertUnsignedCharToFloat(imageUC, image);

	//getting min and max for assignment
	float min, max;
	getMinMaxImage(image, min, max);
	cout << "min of image before filtering: " << min << endl;
	cout << "max of image before filtering: " << max << endl;

  sharpen(image, C);

	getMinMaxImage(image, min, max);
	cout << "min of image after filtering: " << min << endl;
	cout << "max of image after filtering: " << max << endl;

  convertFloatToUnsignedChar(image, imageUC);

  if (!exportFile(imageUC, outputfile)) {
    error(program, "Error writing file");
  }

  cout << "Successfully saved output as " << outputfile << "." << endl;
  return 0;
}

int applyDenoise(int argc, char **argv)
{
  if (argc < 10) {
    error("", "Usage: MoreFiltering -operation Denoise -i <input> -o <degradedoutput> -o2 <blurredoutput> "
      "[-sigma1 s1] [-sigma2 s2] [-f a(veraging)/m(edian)] [-n g(auss)/i(mpulse)]");
  }

  // Store commandline args in stl container
  vector<string> args(argv, argv+argc);

  // Variables to be read from the commandline
	string inputfile, outputfile, outputfileblurred, noise, f, program;
  float s1(2),
    s2(2);

  // Read the arguments.
  getProgramName(args, program);
  getScalarArg(args, "-i", inputfile);
  getScalarArg(args, "-o", outputfile);
  getScalarArg(args, "-o2", outputfileblurred);
  getScalarArg(args, "-sigma1", s1);
  getScalarArg(args, "-sigma2", s2);
  getScalarArg(args, "-f", f);
  getScalarArg(args, "-n", noise);

  cout << "Applying Denoise on " << inputfile << " using sigma1: " << s1 << ", signma2: " << s2 << ", f: " << f << ", n" << noise  << "." <<  endl;

	bool filter(f == "a");

  Image<unsigned char> image;
  if (!importFile(image, inputfile)) {
    error(program, "Error reading image");
  }

  Image<unsigned char> degradedImage(image), degradedBlurredImage(image);

  if(noise == "g") {//gaussian noise
    denoise(image, degradedImage, degradedBlurredImage, s1,  s2, filter, degradeGauss);
  } else {
    denoise(image, degradedImage, degradedBlurredImage, s1,  s2, filter, degradeImpulse);
  }

  cout << "Average Squared Distance (ASD) between original and degraded image = " << averageSquaredDifference(image, degradedImage) << endl;
  cout << "Average Squared Distance (ASD) between original and restored image = " << averageSquaredDifference(image, degradedBlurredImage) << endl;

  if (!exportFile(degradedImage, outputfile)) {
    error(program, "Error writing file");
  }
	cout << "Successfully saved output as " << outputfile << "." << endl;

  if (!exportFile(degradedBlurredImage, outputfileblurred)) {
    error(program, "Error writing file");
  }
  cout << "Successfully saved output2 as " << outputfileblurred << "." << endl;

  return 0;
}

int applyCannyEdgeDection(int argc, char **argv)
{
  if (argc < 9) {
    error("", "Usage: MoreFiltering -operation CannyEdgeDetection -i <input> "
      "-o outputGM, -o2 outputCanny [-s scale]");
  }

  // Store commandline args in stl container
  vector<string> args(argv, argv+argc);

  // Variables to be read from the commandline
  string inputfile, outputfileGM, outputfileCanny, program;
  float scale(1.0f);

  // Read the arguments.
  getProgramName(args, program);
  getScalarArg(args, "-s", scale);
  getScalarArg(args, "-i", inputfile);
  getScalarArg(args, "-o", outputfileGM);
  getScalarArg(args, "-o2", outputfileCanny);

  cout << "Applying Canny Edge Detection on " << inputfile << " using scale: " << scale << "." <<  endl;

  Image<unsigned char> image;
  importFile(image, inputfile);

  cout << "Computing canny edges on scale " << scale << endl;
	Image<float> outCanny;
  cannyEdgeDetector(image, outCanny, scale);

  cout << "Number of edge voxels in canny output: " << nrOfEdgeVoxels(outCanny) << endl;

  Image<unsigned char> outCannyUC(outCanny.dimension());
  convertFloatToUnsignedChar(outCanny, outCannyUC);

  if(!exportFile(outCannyUC, outputfileCanny)) {
    error(program, "Error writing file");
  }

  // Calculate GM image
  Image<float> dxImage, 
    dyImage, 
    gmImage(image.dimension());

  // Compute gradient magnitude and edge direction
  gradientMagnitude(image, dxImage, dyImage, gmImage, scale);

  cerr << "exporting " << outputfileGM << endl;
  Image<unsigned char> outGMuc(gmImage.dimension());
  convertFloatToUnsignedChar(gmImage, outGMuc);
  if(!exportFile(outGMuc, outputfileGM)) {
    error(program, "Error writing file");
  }

  return 0;
}

int main(int argc, char **argv)
{
  if (argc < 3) {
	  error("", "Usage: MoreFiltering -operation Operation [parameters]\n"
      "Only Median, UnsharpMask, Sharpening, Denoise, CannyEdgeDetection are valid operations\n"
	    "Type MoreFiltering -operation Median to get help for the Median operation\n"
	    "Type MoreFiltering -operation UnsharpMask to get help for the UnsharpMask operation\n"
	    "Type MoreFiltering -operation Sharpening to get help for the Sharpening operation\n"
	    "Type MoreFiltering -operation Denoise to get help for the Denoise operation\n"
	    "Type MoreFiltering -operation CannyEdgeDetection to get help for the CannyEdgeDetection operation"
	  );
  }

	string operation(argv[2]);
  if(operation == "Median") {
		return applyMedian(argc, argv);
  } else if(operation == "UnsharpMask") {
		return applyUnsharp(argc, argv);
  } else if(operation == "Sharpening") {
		return applySharpening(argc, argv);
  } else if(operation == "Denoise") {
		return applyDenoise(argc, argv);
  } else if(operation == "CannyEdgeDetection") {
		return applyCannyEdgeDection(argc, argv);
  }

	error("", "Only Median, UnsharpMask, Sharpening, Denoise, CannyEdgeDetection are valid operations");
	return -1;
}